Exchange and correlation in open systems of fluctuating electron number 
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While the exact total energy of a separated open system varies linearly as a function of average 
electron number between adjacent integers, the energy predicted by semi-local density functional 
approximations curves upward and the exact-exchange-only or Hartree-Fock energy downward. As a 
result, semi-local density functionals fail for separated open systems of fluctuating electron number, 
as in stretched molecular ions A^~ and in solid transition metal oxides. We develop an exact-exchange 
theory and an exchange-hole sum rule that explain these failures and we propose a way to correct 
them via a local hybrid functional. 
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Kohn-Sham density functional theory [l], Q (DFT) re- 
places the correlated wavefunction problem by a more 
tractable problem of non-interacting electrons moving in 
self-consistent effective potentials uj(r) (a =f, 1) which 
generate the spin densities n CT (r) of the real (interact- 
ing) system. Exact in principle for the ground-state en- 
ergy and density, Kohn-Sham DFT requires in practice 
an approximation for the exchange-correlation (xc) en- 
ergy functional E xc [n^, rijj. With improving approxima- 
tions, DFT has become the standard method for elec- 
tronic structure calculations in physics and chemistry. 



In terms of the total electron density n 



and 



exchange-correlation energy per electron e xc (r), we write 

E xc = J rfrn(r)e xc (r). (1) 

A ladder 0] of approximations constructs e xc (r) as a func- 
tion of density-dependent ingredients. The first three 
rungs are semi-local (with e xc found from the Kohn- 
Sham orbitals in an infinitesimal neighborhood of r) and 
in some versions non-empirical. The rungs are defined 
by the ingredients: (i) the local spin-density (LSD) ap- 
proximation which uses only n CT (r); (ii) the gener- 
alized gradient approximation (GGA) in the Perdew- 
Burke-Ernzerhof (PBE) Q| version, which adds the gra- 
dients Vn CT (r); (iii) the meta-GGA in the Tao-Perdew- 
Staroverov-Scuseria (TPSS) @ version, which further 
adds the positive orbital kinetic energy densities T CT (r); 
(iv) functionals employing a truly nonlocal ingredient, 
the exact (ex) exchange energy density ne x x , either in 
full [1, 0| or in part 0, HI- The currently used fourth- 
rung functionals are global hybrids (gh) 9], mixtures of 
exact exchange and semi-local (si) approximations 



E: 



gh 



aEZ 



(l-a)E* + E?, (2) 

where the exact- exchange mixing coefficient "a" is a 
global empirical parameter (typically a w 0.2). A global 
hybrid with < a < 1 does not satisfy any universal 
constraints beyond those satisfied by i? xc . 



In many real systems, these existing functionals are 
reasonably accurate for E x and more accurate (due to 
error cancellation) for E xc , with accuracy generally in- 
creasing up the ladder. Yet serious errors occur in nearly- 
separated open systems with fluctuating electron num- 
bers that may not average to integer values, as sum- 
marized below. (a) In the dissociation of heteronu- 
clear diatomics such as NaCl with bond length R, spu- 
rious fractional-charge R — > oo limits are common (e.g., 
Na+ ' 4 • • • Cr - 4 instead of Na° • • • Cl°), as are related 
charge-transfer errors. (In this example, Na- • • CI is a 
closed system of fixed integer electron number, while Na 
and CI are separated open subsystems free to exchange 
electrons with each other), (b) In the dissociation of 
molecular radical cations Aj, the R — > oo limit is cor- 
rectly A+ - 5 • • • A+ - 5 , but the total energy is far below 
that of A- • • A + , with which it should be degenerate. For 
the one-electron molecule H J , this is unambiguously 
a self- interaction error. The errors (a) and (b) are not 
necessarily corrected by functionals that are exact for all 
one-electron densities (c) In the solid state, energy 
competition among electronic configurations in transition 
metal oxides, lanthanides, and actinides can be poorly de- 
scribed [13]. These errors are similar in LSD, PBE and 
TPSS, but are improved by global hybrids [12]. In this 
Report, we derive the generalized x-hole sum rule and 
then explain at the most fundamental level why semi- 
local functionals fail for open systems by showing that 
they violate this rule. We also show that the errors of 
semi-local functionals for E KC can be corrected by a prop- 
erly designed local hybrid functional. 

For closed systems of integer electron number and in- 
teger occupation numbers, the exact (Hartree-Fock-like) 
exchange energy per electron at r is given by 



dr 



, n x (r,r') 



(3) 



where n x (r,r') is the density at r' of the exchange hole 
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FIG. 1: Total energies of the F atom as functions of the aver- 
age electron number N in the PBE approximation (with and 
without correlation) and in Hartree-Fock theory. The dotted 
lines represent the exact-xc result based on the experimental 
ionization potential and electron affinity of the F atom [l5j . 



FIG. 2: Spin-restricted exchange-only dissociation curves for 
the C2 molecule in the symmetric singlet electron configu- 
ration. This configuration is the lowest-energy state in the 
PBEx approximation near the equilibrium internuclear dis- 
tance but not at larger 7?, and not at any R in Hartree-Fock. 



around an electron at r: 
in which 

Mry)=^/ a ^: CT (r')Vw(r) (5) 

a 

is the one-particle density matrix for spin a of the Kohn- 
Sham non-interacting system, tj) acr are orbitals and f aa 
their occupation numbers (f aa = 1 or 0). The density 
for spin a is n a {r) — p a (r, r). 

A little-known aspect of the work of Perdew and 
Zunger [l3[ is their guess that Eqs. ©-(O apply, with 
fractional occupations < f atT < 1, even to an open sub- 
system (with average electron number N = ^ faa) of 
a closed system. We confirm this guess numerically now 
and analytically later. Numerically, we find [l4[ that the 
total Hartree-Fock energy computed self-consistently by 
Eqs. ©-(© for a molecule A 2 at R — > 00 is exactly twice 
that of A +0 - 5 . (The difference between Hartree-Fock and 
exact exchange-only DFT energies is unimportant on the 
scale of the effects we study) . 

Consider the total energy of an F atom, treated 
as an open system, as a function of average electron 
number N. We vary N by changing the population 
of the highest partly-occupied 2p-orbital and compute 
the e nerg y self-consistently using the GAUSSIAN pro- 
gram [16|. Fig. [T] shows that the total energy with exact 
exchange-correlation varies linearly between adjacent in- 
tegers (3, [13, but the PBEx and PBExc energies curve 
upward strongly, while the HF energy curves downward. 
Note also the accuracy of PBExc energy differences for 
integer values of N around 9 and their inaccuracy for 



non-integer values. LSD [l7| . TPSS [HI, and other func- 
tionals [X 81 ] behave similarly. 

In Fig. [TJ the Hartree-Fock approximation (exact ex- 
change without correlation) shows substantial midpoint 
error for non-integer electron numbers. However, this er- 
ror is positive and thus not very harmful since energy 
minimization forces integer electron numbers onto sepa- 
rated open subsystems, e.g., a symmetry-broken F- ■ • F + 
as the dissociation limit of But F+ 5 • ■ ■ F+ 5 is a 
harmfully too-deep minimum for the semi-local density 
functionals, since their midpoint error is negative. 

The failure of semi-local exchange can occur even in 
open systems with integer average electron number. In 
semi-local DFT, the correct dissociation limit for neutral 
homonuclcar diatomics built up from open-shell atoms 
(e.g., H2) is achieved by spin symmetry breaking. If, 
however, spin symmetry is imposed, then the semi-local 
exchange energy in the R — > 00 limit is much more nega- 
tive than the exact exchange energy [l9j . Typically, the 
separated atoms have half-integer numbers of electrons 
of each spin (e.g., N^=N^ = ^ on each spin-unpolarized H 
atom) . But the C2 molecule in the unusual singlet config- 
uration KK(<jg2s) 2 (cr u 2s) 2 (ir u 2p x ) 2 (Tr u 2py) 2 dissociates 
to two neutral C atoms, each with iVj=iVj=3 but frac- 
tional occupations f 2p ^=f 2pA =f 2py ^=f 2pyl = l. Fig. [5] 
shows that, again, the spin-restricted semi-local exchange 
energy in the R — * 00 limit is far more negative than that 
of spin-restricted Hartree-Fock theory. 

In a transition-state complex of a chemical reaction, 
residual fluctuation of electrons among its weakly-bonded 
fragments raises the total energy via increased Coulomb 
repulsion, but semi-local exchange approximations miss 
this effect and predict reaction barriers that are too low. 

When we apply single-configuration Hartree-Fock the- 
ory to a closed system, fractional occupation numbers on 
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a separated open subsystem can only represent fluctua- 
tion of electrons among separated subsystems. That is 
the case we address here. In a different case, fractional 
occupation on a closed system can represent fluctuation 
of electrons among degenerate orbitals of that system. 

We will now prove that Eqs. ©-(O apply to an open 
system of fluctuating electron number and then show how 
they imply the behavior found in Figs. Q] and O Let S 
be a fully separated open system and R its reservoir, i.e., 
the only other open system with which S is free to ex- 
change electrons. Consider now the closed system S+R 
of integer electron number, which we will describe by a 
single Slater determinant. If the wavefunction of S+R 
were fully correlated, the ensemble describing S as a sub- 
system of S+R would also be [T3|, but since the former is 
only a single energy-minimized Slater determinant, the 
latter consists of Slater determinants and their proba- 
bilities that need not be energy-minimized. When we 
form S+R, each orthonormal orbital on S remains un- 
changed if there is no orbital of the same energy from R, 
or hybridizes (forms a linear combination) with orbitals 
of the same energy from R, resulting in a band whose 
width tends to zero in the infinite-separation limit (i.e., 
a group of nearly degenerate molecular orbitals) . If this 
band is filled, hybridization is just a unitary transforma- 
tion of the occupied orbitals. The Hartree-Fock exchange 
energy of S+R, as predicted by Eqs. (JS])-©, is a sum of 
a pure S term (arising when r and r' are on S) and a pure 
R term (arising when r and r' are on R), because when r 
is in S and r' in R the Coulomb interaction l/|r — r'| van- 
ishes. Thus, the S term is given by Eqs. ©-(O in which 
ifi aa - are the orthonormal orbitals of S with f aa = 1 if 
the corresponding band is filled and < j aa < 1 if it is 
partly filled (holds more nearly-degenerate orbitals than 
electrons), as it can be if the Fermi level lies within the 
band. This concludes the proof. 

Now we can use the ensemble describing S to explain 
why the Hartree-Fock energy as a function of N curves 
downward. Let £' l HF (S; S+R) be the Hartree-Fock energy 
of the TVi-electron pure state i of system S, evaluated from 
orbitals for S formed by truncating and renormalizing the 
ground-state molecular orbitals of S+R, and let pi be the 
probability to find state i in the ground state of S+R, 
with Y^iPi = 1 and J2iPi N i = N - Then # HF (S; S+R) = 
EiPiE? F (S;S+R) > EiPi#? P (S;S). A familiar and 
instructive example 20] is the H atom in a spin-restricted 
stretched H2 molecule, where the states i are the neutral 
atoms of each spin (pi — 0.25, Ni = 1), cation (pi = 
0.25, Ni = 0), and anion (pi = 0.25, N = 2). Another 
example [2l[ is Fig. [1] where the states i are those for 
Ni = J and J — 1 electrons, with J an integer. 

The reason why semi-local functionals predict too- 
negative energies for systems with fractional occupations 
also becomes clear. By the orthonormality of the orbitals 
of a subsystem, Eqs. ((4]) and ([5|) imply the sum rule 

/dr'n x (r,rO = -EE /i ^ M ' ( 6 ) 



where n aa (r) — f a a\ipaa(^)\ 2 ■ Eq. © fo r noninteger 
f aa was presented without proof in Ref. [lj. Adding and 
subtracting — 1 = — J2aa n aa (Y)/n{Y) to Eq. (J6j) , we find 

/ dr> n x (r, r') = -1 + V / QCT (1 - f aty ) ^ a °^ ? . (7) 
J ^ n(r) 

The close similarity between the sum rule integral [left- 
hand side of Eq. 0] and the exchange-energy-per- 
electron integral [right-hand side of Eq. ([3])] for integer 
f aa was pointed out by Gunnarsson and Lundqvist [22|, 
and this similarity persists for noninteger values of /„„. 
When all the occupation numbers are 1 or 0, the right- 
hand side of Eq. ([7]) becomes — 1, which is also the sum 
rule implicitly assumed by LSD, PBE, and TPSS [H[. 
But, when some occupation numbers are between 1 and 
0, the right-hand side of Eq. ([7]) will fall between —1 and 
0, meaning that part of the exact exchange hole around 
an electron in an open system is located in its distant 
reservoir. In this case a semi-local exchange approxima- 
tion e x will tend to be more negative than the exact ex- 
change energy per electron e x x (r), as illustrated in Figs. [I] 
andO The total or system-averaged £" x also violates the 
sum rule on the system-averaged [24| x-hole 

i f dr'drn(rK(r,r') = -l + ^^/ QCT (l-/ QCT ). (8) 

J eta 

The Hartree-Fock mean square fluctuation of electron 
number, ^ Q f aff (l — faa), only vanishes when all oc- 
cupations are or 1. 

A sum rule for the exact xc-hole density rt xc (r,r') is 
also known [l7|. Its integral equals —1 only when the 
electron number does not fluctuate and otherwise falls 
between —1 and 0. Ref. presents a coupling-constant 
integration for E xc and n xc . But the integrand for E xc at 
zero coupling strength is not really the exact exchange- 
only energy because of an exact-degeneracy static cor- 
relation. When the electron number on the infinitely- 
separated open system S fluctuates at the Hartree-Fock 
level, occupied and unoccupied orbitals (with the same 
spin) of closed system S+R are degenerate. Degenerate 
perturbation theory is needed to find the correlation en- 
ergy, which is of the same order as the exchange energy 
even in the weak-coupling or high-density limits. Exact- 
degeneracy correlation and normal correlation shift the 
downward-curved Hartree-Fock energy of Fig. Q] into the 
straight-line correlated exact energy. Note also from 
Fig. [T] that semi-local approximations for E xc overesti- 
mate the strength of exact-degeneracy correlation (which 
they introduce via E^r). 

Semi-local functionals are often combined with a Hub- 
bard U parameter (DFT+Z7). A simple case occurs when 
only one localized orbital has non-integer occupation /, 
and the method adds to the semi-local energy a positive 
term Uf(l — /). The close connection between DFT+?7 
and self- interaction correction has been argued [25[ . We 
note that U does not represent "strong correlation" (as 
sometimes asserted), because the U needed to reach the 
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Hartree-Fock energy is greater than that needed to reach 
the exact correlated energy. U favors the less fluctuating 
configuration by penalizing the more fluctuating one. 

Our Fig. [1] and our analysis explaining it show 
that some region-dependent fraction (between and 
1) of exact exchange is needed to correct semi- local 
exchange-correlation approximations. Such a mixing of 
the downward-curving exact exchange with the upward- 
curving semi-local exchange and semi-local correlation 
can produce the needed straight line. This motivates a 
local hybrid (lh) functional 

e> h c (r) = + [1 - a(r)](4 - ef ) + e* (9) 

where < a(r) < 1 and sl=TPSS. Eq. ([9]) was presented 
in Ref . J2H without a form for a(r). Forms were proposed 
in Ref . y and in Ref . |7| (where the term "local hybrid" was 
coined), but did not achieve useful accuracy for equilib- 
rium properties Q. The choice a(r) = 1 satisfies nearly 
all universal constraints but misses the delicate and help- 
ful error cancellation between semi-local exchange and 
semi-local correlation that typically occurs (because the 
xc-hole is deeper and more short-ranged than the x-hole) 
in normal regions where the density is not too high, too 
strongly- varying, too one-electron-like, or too fluctuating 
in a spin-polarized region at the Hartree-Fock level (i.e., 
e xV e x ^ So we will take o(r) to be small or in a 
normal region, and to tend toward 1 to the extent that 
any condition of normality is violated. Then the second 
and third terms on the right of Eq. ([9]) represent [l9| 
near-degeneracy static and dynamic correlation, respec- 
tively. The dominance of exact exchange [a(r) — » 1] in 
the high-density limit means that Eq. © has 100% exact 
exchange plus fully nonlocal correlation. 



A natural generalization of the global hybrid of Eq. (J2J, 
Eq. © can satisfy [3|, |27| many more exact constraints 
while achie ving greater accuracy. We have developed and 
are testing [27| such a local hybrid hyper-GGA, and the 
results will be reported later. One or more universal em- 
pirical parameters are needed, as in Eq. since the 
universal constraints are already satisfied by a(r) = 1. 
As usual, symmetry must be allowed to break. Spin- 
symmetry breaking [281 ], like real correlation, lowers the 
energy by suppressing Hartree-Fock-level fluctuation of 
electron number. Semi-local functionals mimic this sup- 
pression. Where real correlation cannot fully suppress 
fluctuation, as in open systems of non-integer average 
electron number, semi-local functionals overcorrelate and 
need a large positive correction, i.e., a large a(r). 

The Perdew-Zunger self-interaction correction [l3[ to 
semi-local functionals works in much the same way to 
raise the energy of a system with fractional occupa- 
tion [llj], satisfying the sum rule of Eq. (0. However, it 
loses the error cancellation between semi-local exchange 
and semi-local correlation in normal regions, and so is in- 
accurate for molecules near equilibrium geometries [29j. 

In summary, striking and diverse failures of semi-local 
functionals arise because they assign too low an en- 
ergy to configurations where the electron number in a 
spin-polarized region fluctuates too strongly (i.e., where 
e xV e x ^ 1) a * the Hartree-Fock level. These errors 
should be corrected by local mixing of a variable frac- 
tion of exact exchange. 
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